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Abstract 

A deterministic model of an age-structured population with genetics 
analogous to the discrete time Penna model (TJ [2] of genetic evolution is 
constructed on the basis of the Lotka-Volterra scheme. It is shown that if, 
as in the Penna model, genetic information is represented by the fraction of 
defective genes in the population, the population numbers for each specific 
individual's age are represented by exactly the same functions of age in 
both models. This gives us a new possibility to consider multi-species 
evolution without using detailed microscopic Penna model. 

We discuss a particular case of the predator-prey system representing 
an ecosystem consisting of a limited amount of energy resources consumed 
by the age-structured species living in this ecosystem. Then, the increase 
in number of the individuals in the population under consideration de- 
pends on the available energy resources, the shape of the distribution 
function of defective genes in the population and the fertility age. We 
show that these parameters determine the trend toward equilibrium of 
the whole ecosystem. 



1 Introduction 

Thousands of papers have been published on the Lotka-Volterra equations [H [4] 
describing population growth, competition or speciation. In real populations 
the reproduction rate of individuals depends on their age and therefore it is 
necessary to include age structure into these equations. One example of how 
this can be achieved can be found in [5]. It is also possible to introduce a time 
delay between cause and effect (see, e.g. [6]). However, the majority of the 
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Figure 1: Monte Carlo simulation results of the Penna model of genetic evolu- 
tion after 5 x 10 7 MC steps. The upper panel shows the age-specific fraction 
of defective genes in three haploid populations with different values of the in- 
dividual's fertility age of = 15, 30, 40. The remaining parameters of the Penna 
model are the same, T = 1, M = 1, bit-string length n — 100. The resulting 
limit age of individuals in these populations ao = 29,45,55 is the age where 
the fraction of defective genes is equal to 1. The lower panel shows the age 
distribution in these populations. 



Lotka-Volterra equations do not usually include genetic information and the 
question arises how to include it directly into the Lotka-Volterra equations. 

The age-specific equations for population growth seem to be a good candi- 
date to represent genetic information because the age structure introduces some 
analogy to the Penna model [TJ [2] of genetic evolution. It is a model of genetic 
evolution although all details concerning the genes are skipped except of the 
state of their functionality - is the gene under consideration correct or it is mu- 
tated. This simple model has turned out to be very successful in interpreting the 
demographers data of real populations even as complex as human populations 

USEE]. 

In the original asexual version of the Penna model p] , the population under 
consideration consists of individuals represented by genomes defined as a string 
of n bits. The bits represent states of genes where denotes its functional allele 
and 1 its bad allele. It is assumed that if an individual possesses T bad alleles 



2 




5 10 15 20 25 30 5 10 15 20 25 30 



a a 

Figure 2: The dependence of the population number on the individual's age in 
the case of the equilibrium Penna model of the bit-string evolution (T = l,a_p = 
15, M = 1) and the age-structured Lotka-Volterra equations (the left panel of 
the figure). The right panel shows the distribution of defective genes in the 
equilibrium Penna model for the whole population and for individuals of age 
> Of. The solid lines represent an approximation of the data e a (fraction of 
defective genes in the population) which is power-like for the ages a < ap and 
exponential-like for a > aF. 

switched on, it dies. In the model all genes are switched on chronologically - 
each bit corresponds to one "year" - and maximum life span of an individual 
is ao "years". After reaching the fertility age af an individual gives birth to 
B offsprings whose genomes are mutated versions of the parental genome. The 
mutation rate M is constant. After the mutation the affected gene is represented 
by a bit with a value opposite to the value before the mutation. The results of the 
Monte Carlo simulations of three different haploid populations for the Penna 
model in the case when back mutations from 1 to zero are not allowed have 
been shown in Fig. [I] In the figure, there have been plotted fraction of defective 
chronological genes in each population and age distribution of individuals. 

In the diploid version of the Penna model the individual's genome is repre- 
sented by two bitstrings and then each locus possesses two alleles. The diploid 
individuals can reproduce sexually [2]. Both the haploid model of the equilib- 
rium population and its diploid version are uniquely described by the fraction 
of defective genes in the population specific for each individual's age. A short 
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review of Monte Carlo simulation results for the Penna model with some addi- 
tional details, like the presence of the housekeeping genes or the recombination 
frequency, can be found in [S]. 

Although the age distribution curves obtained in the Penna model coincide 
very well with demographers data for real populations, the model is not so "in- 
terdisciplinary" as the Lotka-Volterra population model. The reason could be 
that it is very difficult to obtain analytical results for the general Penna model 
and therefore one has to use the Monte Carlo method. However Monte Carlo 
simulations of the Penna model need large populations and the simulation time 
is very long. Hence, the typical multi-species problem exceeds the computing 
capability of a single PC-computer. Another problem is how to avoid the corre- 
lations arising from parallel computing if one tries to distribute the simulations 
to many processors. On the other side hand, it is relatively easy to solve numer- 
ically even a large set of the differential equations describing the Lotka-Volterra 
populations. Below we show how to include the fraction of defective genes into 
the equation for population growth so that the age distribution curves for these 
two models coincide. 

In the following sections we discuss the behavior of the Lotka-Volterra ecosys- 
tem in which the age-specific species has been determined by the form of the 
fraction of defective genes. 

2 Population growth of a single species 

In this section we restrict ourselves to the haploid version of the Penna model 
but the results could be generalized to diploid populations. All we need from the 
Penna model is the distribution of defective genes in an equilibrium population 
and the corresponding age distribution. Such data could be taken from a real 
population as well. In the case of the Penna model we have performed a series 
of simulations of genetic evolution of the bitstring populations for different pa- 
rameters like the fertility age a_F, genome length and the number of offsprings 
born each year. We considered the simplified case when the parameter T = 1. 
The Verhulst factor was used to control the population size. 

In the right panel of the Fig. [2] in the semi- log scale there have been plotted 
separately the fractions e a of the age-specific defective genes in the whole equi- 
librium population and in the part of it which is consisting of the individuals 
with the age a > ap. Fraction e a of the age-specific defective genes in the whole 
equilibrium population and in the part consisting of individuals with a > ap 
are plotted in Fig. [5] on semilog axes. In the latter case all individuals should 
posses good genes specific for a < ap (e a — 0) because T = 1. Otherwise they 
should have died. 

Our deterministic model for population growth is constructed in such a way 
that the values of the fractions e a for each age a = 0, 1, . . . , ac are inserted into 
the following set of the age-specific differential equations describing the change 
in number of individuals of age a: 
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= a a _iV(t)(l - e a _i)JV _i(t) - j a N a (t), (1) 
where JV a is the number of individuals specific for the age a = 1, . . . , oq 

iv = x; *«(*)■ ( 2 ) 

a>aF 

and 1/ represents the Verhulst saturation factor 

V (t) = l-iY,N a (t). (3) 

a=l 

with 57 being the saturation level. The first term is the graduation term from 
the age a — 1 to the age a whereas the other part of the equation describes the 
graduation to the age a+1 and the individual's mortality at age a. These terms 
are controlled by the age-specific rate coefficients a a and j a , respectively. If we 
restrict the above to genetic evolution only, as in the case of the Penna model, 
we would have to choose the values j a = 1 for all ages a, the reproduction rate 
cto should be chosen as in the computer simulations of the Penna model (e.g. 
we set it to 1.1) but a a = 1 for all ages a > 0. The Verhulst factor (the same as 
in the Penna model) ensures that the solution of the equations Eq. (TTJ) saturates 
at long times and then the age profile in the population coincides with the one 
from the equilibrium Penna population. An example of this can be observed in 
the left part of Fig. [21 where an analytical approximation of the simulation data 
{£ a }i D has been applied. 

Equations Eq.{T]) seem to possess many parameters. However, the simula- 
tions of the simple haploid version of the Penna model for different values of clf 
suggest that in the case of T = 1 the fraction of the age-specific defective genes 
in the population which can be activated, i.e. after which the individuals die, 
consists of two parts (Fig. [2]), one part relates to individuals of age a > of and 
another one relates to individuals of age a < clf- In this simple case we can 
approximate them as follows: 



for a < clf and 



Ca d (4) 



Ae ba (5) 



for a > cif, with constants A, b, C, d. From the normalization condition, 



Ea D = 1, we obtain 



A = e- baD . (6) 

Then, the number of parameters describing the fractions of defective genes drops 
down. It can be even smaller if we continue the left branch of the function s a 
to a = a,F + 1 and merge it with the right branch at a = ap + 1 in such a way 
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that individuals of age a > did not contributed into the left branch also at 
a = ap + 1, i.e. 

C(l + a F ) d - Ae 1+aF = Ae 1+aF . (7) 
Then we could estimate the value of the parameter C: 

C = 2Ae b{1+aF) /(l + a F ) d . (8) 

In this way we have got an approximate distribution of defective genes controlled 
by four parameters, b, d, ap and ap>- In the equations Eq.(TT]) there is also the 
parameter ao, but it was kept constant and its value was the same as in the 
computer simulations of the Penna model. 

The small number of parameters in the deterministic model makes it pos- 
sible to consider multi-species evolution where the species differ in the values 
of of, do, b, d and ao. In this case the numerical solutions do not require a 
high computational effort as in the Penna model. We should add that if the 
distribution function of the defective genes is an unknown function of a and 
only the empirical values e a are available, then these values could be used in the 
same way as in the discussed example, i.e., we do not need the analytical form 
of the distribution function for defective genes. We have considered only the 
case of the genetic evolution in equations Eq.|T]) but there is a rich possibility 
to use other values of the coefficients a a and 7 a if one wants to describe real 
populations. 



3 Predator-prey system 

We consider an ecosystem in which a limited amount of the self-regenerating 
energy resources, S(t), are consumed by the age-structured species living in this 
ecosystem. In this case, the increase in number of the individuals in the species 
depends on the available energy resources necessary for life processes. We will 
restrict ourselves to one species only. The species will carry genetic information 
represented by the fraction of defective genes in the population, the fertility age 
and the longevity ap, . We will show that the change of these genetic parameters 
can change the trend toward equilibrium of the whole ecosystem. 

In the case of one species the ecosystem can be described by the following 
predator-prey set of equations 



^ = as S(t){l-^)- ls N(t)Sit) 



dN a (t) 
dt 



= a a -iS(t)(l-e a -i)N a -i(t)-N a (t), 

(9) 



where a — 1, . . . , ap>, Nq is the same as in Eq.Q, N(t) is the total number of 
individuals 



G 



a=l 



(10) 



S(t) represents energy units consumed by the species, fls represents the satura- 
tion level for the energy resources, the coefficient as is the regeneration rate of 
the ecosystem resources, 7s is the damping coefficient due to the species using 
the energy resources. In order to follow the Penna model of genetic evolution 
the species evolves according to the equations Eq.fl]) (*y a = 1 for all values of 
a) but the species' growth is controlled by the amount S(t) of the available 
energy resources instead of the Verhulst factor, V(t). However, a Verhulst term 
is assumed for the self-regeneration of the ecosystem with respect to the energy 
resources. The energy resources represent prey and the species plays the role 
of predator. One could easily generalize the above set of equations to include 
many species. 

We could expect that equations Eq. ([9]) have solutions which change in time in 
a similar way as the microscopic Penna model evolution in a limited ecosystem. 
We have recently studied an analogous predator - prey problem in [lOj , where a 
variable surrounding and its effect on the Penna model was considered. One of 
the findings of the computer simulations was that small changes in the inherited 
genetic information could lead to spontaneous bursts of evolutionary activity. 
The predator-prey dynamics with genetics could be as complex as discussed by 
Ray et al. pj] who showed that the system passes from the oscillatory solution 
of the Lotka-Volterra equations into a steady-state regime, which exhibits some 
features of self-organized criticality (SOC). 

In our model, we can always expect saturation of the energy resources if 
there is no species (N(t)=0) because the Verhulst term is present in Eq.©. 
However, if N(t) > then depending on the shape of the distribution function 
of defective genes the solution of Eq. © saturates as in the example in the left 
part of Fig. [3] or it oscillates as in the right part of the figure. This means that 
an ecosystem with a few species can exhibit very complex evolution in which 
the solutions representing some species will spiral towards an equilibrium fixed 
point and some of them will try to converge to the limit cycle. The smaller 
the value of of is, the larger the oscillations are, and the species may become 
extinct. 

It is important that not only genetic parameters can change the type of the 
asymptotic solutions of Eq.©. Consider a species similar to the one represented 
in the left hand part of, Fig. [3] for which the Lotka-Volterra solutions spiral to 
a fixed point. Let's increase the value of the fertility age from clf to a' F > ap 
after some period of time but let the genetic information represented by the 
values {6a}i D remain the same as before. For example, according to some social 
regulations since a specific time moment the individuals start to reproduce at 
older ages a > a pi = 20. The consequences of this shift in the reproduction age 
is that the solutions of Eq. (J9j) change qualitatively from the trend approaching 
the fixed point as in the left part of Fig. [3] to oscillations as in the case of the 
species from the right hand part of Fig. [3j It could happen that such a transition 
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Figure 3: Both panels show the variation of the energy resources S(t) and the 
total number of individuals N(t) of the age structured species starting from the 
same initial conditions. The difference is that on the left clf — 15 and on the 
right a_F = 20. The remaining parameters are the same, a_o = 29 and ao = 1.1. 



between these two types of solutions of the Lotka-Volterra equations could cause 
the species to become extinct. 

In [T2] it has been shown that the Lotka-Volterra systems have a self- 
regulatory character and there exist threshold values for the fraction of destroyed 
population above which the system returns to its previous state. We observe a 
similar behavior in our model. In particular, in Fig. @] an effect of three such 
catastrophes on two different species represented by two different types of solu- 
tions of the Lotka-Volterra equations is shown. It is interesting that for some 
solutions the introduction of a perturbation decreases the amplitude of their 
oscillations. This could be observed in Fig. HI The trace of this catastrophe in 
the (S(i), N(t)) representation has been shown in Fig.[5j 

Thus, we could expect that in the ecosystem consisting of many species a 
rapid change like the extinction of a species or a sudden increase in their number 
can cause discontinuous changes of the fluctuations of the population numbers 
of the existing species, similarly as in Fig. [5] or Fig. [5] (discussed below). The 
age-specific structure of the competing populations seems to be an important 
factor in their survival. 

It is easy to generalize the Lotka-Volterra equations Eq.([9]) to the case of 
many species competing for the resources. If index k runs through M different 
species, the age-structured growth equations should be changed to the following: 

rlN k (t\ 

— £± = a*_ x S(t)(l - e*_i)JVti(*) - N a(t), (11) 
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Figure 4: The effect of three catastrophes, where 30% of the population under 
consideration has been eliminated at two time moments t — 300 and 400 and 
10% of the population at t = 500. On the left there is a species with ap = 15 
and od = 29 and on the right the species has op = 20 and an = 29. In both 
cases the catastrophe at t = 500 decreased the population number fluctuations. 
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Figure 5: The effect of catastrophe from Fig.|4]on the variation of S(t) and N(t) 
in the population with oj? = 20 and ao = 29. 
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Figure 6: Dependence of population number N(t) on time t in the case of two 
species (a) and (b) competing for the same energy resources. In both figures, 
the species (a) has the same fertility age ap = 15 but an = 23 on the left and 
an = 24 on the right. The species (b) (ap = 20 and an = 29) is the same in 
both cases. The same initial condition have been applied. 

where fc = 1,2,..., M. Even in the case of two species M = 2 the ecosystem 
under consideration can exhibit complex behaviour. For example, Fig. [S] shows 
the time dependence of the population number N(t) for two competing species 
(a) and (b) in two cases which differ for species (a) by one parameter only, 
an = 23 (left) and ap, — 24 (right). The shift in the value of an decides which 
species may become extinct. 

4 Conclusions 

It has been shown how to construct age-structured Lotka-Volterra equations in 
order to describe some features of genetic evolution. This could be helpful in 
modeling real populations where many parameters are typically used. In the 
considered model all species were competing for the same energy resources. 
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